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ABSTRACT 

We use data from the first epoch of observations with the Infrared Array 
Camera (IRAC) on the Spitzer Space Telescope for the Great Observatories Ori- 
gins Deep Survey (GOODS) to detect and study a collection of Lyman-break 
galaxies at 2 ~ 6 to 5 in the Hubble Ultradeep Field (HUDF), six of which have 
spectroscopic confirmation. At these redshifts, IRAC samples rest-frame optical 
light in the range 0.5 to 0.8/im, where the effects of dust extinction are smaller 
and the sensitivity to light from evolved stars is greater than at shorter, rest- 
frame ultraviolet wavelengths observable from the ground or with the Hubble 
Space Telescope. As such, it provides useful constraints on the ages and masses 
of these galaxies' stellar populations. We find that the spectral energy distribu- 
tions for many of these galaxies are best fitted by models of stellar populations 
with masses of a few xlO^^ Mq, and with ages of a few hundred million years, 
values quite similar to those derived for typical Lyman break galaxies at 2 ~ 3. 
When the universe was only 1 Gyr old, some galaxies had already formed a mass 
of stars approaching that of the present-day Milky Way, and that they started 
forming those stars at z > 7, and in some cases much earlier. We find that the 
lower limits to the space density for galaxies in this mass range are consistent 
with predictions from recent hydrodynamic simulations of structure formation in 
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a ACDM universe. All objects in our samples are consistent with having solar 
metallicity, suggesting that they might have already been significantly polluted 
by metals and thus are not comprised of "first stars". The values for dust red- 
dening derived from the model fitting are low or zero, and we find that some of 
the galaxies have rest-frame UV colors that arc even bluer than those predicted 
by the stellar population models to which we compare them. These colors might 
be attributed to the presence of very massive stars (> IOOMq), or by weaker 
intergalactic HI absorption than what is commonly assumed. 

Subject headings: cosmology: observations — galaxies: evolution — galaxies: 
luminosity function, mass function — infrared: galaxies 

1. Introduction 

A striking result from studies of galaxy formation and evolution is that galaxies as 
massive as ~ existed at z 1 or even out to z ^ 2 (e.g., Fontana et al. 2004, 

Cimatti et al. 2004, Glazebrook et al. 2004, Daddi et al. 2005), and galaxies with masses 
~ IO^^Mq are already rather commonly seen at z ^ 3 (e.g., Papovich, Dickinson & Ferguson 
2001; Shapley et al. 2001). Many such massive galaxies apparently have well-evolved stellar 
populations, i.e., they must have formed their stars well before the epoch at which they 
are observed. For example, about 20% of the LBGs in the z ^ 3 sample of Shapley et al. 
(2001) have stellar masses > lO^^M© and inferred ages > IGyr, implying formation redshifts 
oi Zf > 5. Infrared selected samples have identified candidates for more massive, evolved 
galaxies at, z ^ 2-3 (Franx et al. 2003; Yan et al. 2004; Labbe et al. 2005), with inferred 
stellar masses that can exceed IO^^Mq. The dominant luminosity components for some of 
these objects are best explained by an old stellar population with ages of 1.5-2.5 Gyr, which, 
taken at face value, suggests that they formed no later than Zf — 5, and possibly as early as 
Zf = 15-20. 

If the inferences from these studies are correct, we should see massive galaxies at z ^ 5 
and beyond. One of the central science drivers of the Great Observatories Origins Deep 
Survey (GOODS) Spitzer Legacy Science Program is to constrain the stellar masses of L* 
galaxies out to 2; ~ 5 with the rest-frame ncar-lR observations made by the Infrared Array 
Camera (IRAC; Fazio et al. 2004). In fact, the performance of this instrument in its 3.6 
and 4.5//m channels exceeds the pre-launch expectations because the in-fiight point spread 
function (PSF) exceeds the telescope specifications. This enables the detection of galaxies as 
distant as 2; ~ 6. At these redshifts, the observed optical and near-infrared data sample rest- 
frame ultraviolet light, which is sensitive to young star formation and dust absorption, but 
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which gives relatively little information about the age or mass of a galaxy's stellar population. 
Spitzer and IRAC offer access to the rest-frame optical light which is less subject to dust 
extinction and more sensitive to the longer-lived stars that dominate the stellar mass. Egami 
et al. (2005) have used IRAC to detect a galaxy at 2; ~ 6.7 Icnscd by a foreground cluster, 
while Eyles et al. (2005) have studied two galaxies at z ^ 6 that are detected in the GOODS 
IRAC data. In both papers, the spectral energy distributions were used to estimate the 
stellar masses of these galaxies, and to demonstrate the presence of stars formed several tens 
to several hundred million years earlier. 

In this paper, wc study the rcst-framc UV to optical photometry of 2; ~ 6 and 5 (hereafter 
"high-z") galaxies selected from the Hubble Ultra Deep Field (HUDF^; PI. S. Beckwith), 
using the HST ACS/NICMOS data and the first epoch of GOODS IRAC observations of 
the Chandra Deep Field South (CDF-S). While the GOODS IRAC data cover a much wider 
area, an important advantage of analyzing galaxies in the HUDF is the availability of very 
deep near-infrared imaging, in two bands, from HST/NICMOS (Thompson et al. 2005) for 
a significant number of objects. Without such data, there is rather little constraint on the 
rest-frame ultraviolet properties of galaxies at z > 5 compared to what other studies have 
had available when analyzing Lyman-break galaxies at lower redshifts. The UV luminosity 
of a galaxy is an important constraint on its rate of on-going star formation, while the UV 
spectral slope is sensitive to dust reddening. When fitting population synthesis models to 
galaxy photometry, it is important to have constraints on both of these quantities in order 
to derive reliable star formation rates, extinction, stellar population ages, and mass-to-light 
ratios. At 2; ~ 6, flux in the ACS ^775 band, and to a lesser extent also the ZgsQ-band 
(depending on the galaxy redshift), is suppressed by the stochastic effects of the Lyman a 
forest. The ACS bands alone are insufficient for reliably measuring the UV luminosity of 
z ^ 6 galaxies, and the UV slope cannot be constrained without redder bandpasses. At 
2; ~ 5, the 2:850 band gives a reliable luminosity at ~ 1500A, but the spectral slope is again 
unconstrained without near- infrared data. These high-rcdshift galaxies are extremely faint, 
and only the brightest are detected even in the very deep near-infrared data available for the 
GOODS project from ISAAC on the VLT, and even then at fairly low signal-to-noise ratios. 
The NICMOS HUDF data are much deeper and detect all galaxies from our sample that fall 
within their field of view, providing two more photometric bands that cleanly constrain the 
UV luminosity and spectral slope. 

The paper is organized as follows. The high- 2; galaxy samples are briefiy described in 
§2, and the IRAC counterparts of these galaxies are discussed in §3. The stellar population 
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synthesis models that we use to fit the observed spectral energy distributions (SEDs) of these 
galaxies are described in §4. The constraints on the stellar populations at z ~ 6 and 5 are 
given in §5 and 6, respectively. We briefly compare our results with recent hydrodynamic 
simulations in §7, and conclude the paper with a summary in §8. We denote the F435W, 
F606W, F775W, and F850LP bands of ACS as -B435, Vqqq, ijy^, and zs5o, respectively, and 
the FllOW and F160W bands of NICMOS as Jno and ifieo, respectively. In one occasion 
we also use the X^-band photometry obtained by the ISAAC instrument at the VLT. All 
magnitudes are in the AB system. Throughout this paper, we adopt the following cosmolog- 
ical parameters based on the Wilkinson Microwave Anisotropy Probe result (WMAP) from 
Spergel et al. (2003): Qm = 0.27, Qa = 0.73, and Ho = 71 kms-^ Mpc"^ {h = 0.71). 

2. Samples of Galeixy Candidates at z ^ 6 and 5 in the HUDF 

The sample of z ~ 6 galaxy candidates used in this paper is from Yan & Windhorst 
(2004; hereafter YW04), which consists of 108 objects to Zs5o < 30.0 mag. These candidates 
were selected as iyys-band dropouts using the criteria of (1) (^775 — 2^850) > 1-3 mag and (2) 
non-detection in both the S435 and the Veoe bands {S/N < 2). These criteria are similar to 
those used in Dickinson et al. (2004). The targeted redshift window is 5.5 ^ z 6.5. There 
are 6 multiple systems among these candidates (see Table 1 of YW04), where a multiple 
system is defined as a group whose members are within 1" of each other. Such a multiple 
system cannot be resolved by IRAC, and will show up as a single IRAC source if detected. 

The sample of 2; ~ 5 galaxy candidates is the Vgoe-band dropout sample of Yan et al. 
(in preparation). Briefly, the selection criteria are (1) 



and (2) non-detection in the ^435 band {S/N < 2). These criteria are similar to those of 
Giavalisco et al. (2004b), but are fine-tuned to better suit the high S/N HUDF data, which 
have a tighter low-redshift galaxy locus in the color space. The targeted redshift window is 
4.5 < 2; < 5.5. This sample consists of ~ 550 candidates, among which 95 objects form 43 
multiple systems. 




(1^606 - ^775) > 0.4 + 2.67 X (i775 - Zs5o), 



if (^775 - Zs5o) < 0.3 

if 0.3 < (i775 - Z85o) < 0.6 

if (^775 - Zs5o) > 0.6 
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3. IRAC Counterpcirts to High-z Galstxy Candidates 

The IRAC data used in this paper are the mosaics of the first epoch of GOODS observa- 
tions of the CDF-S ^. These data are discussed in detail in Dickinson et al. (in preparation), 
and are also briefly described in Yan et al. (2004). To summarize, the middle one-third of 
the GOODS GDF-S field, which includes the HUDF, has been observed in all four IRAG 
channels with a nominal exposure time of ~ 23.18 hours per pixel. The final drizzle-combined 
mosaics have a pixel scale of 0.6", or approximately half of the native IRAG pixel size. These 
IRAG mosaics are registered to the same astrometric grid as the GOODS AGS images of the 
same area, which is also the astrometric frame of the ACS HUDF. Sources were detected in 
a weighted sum of the 3.6;um and 4.5/im images. The IRAC photometric catalogs used here 
are somewhat different from those used in Yan et al. (2004) in the sense that the source 
extraction parameters for the current catalogs are better tuned for deblending objects. 

3.1. Matching high-z objects in the IRAC images 

As the FWHM of the PSF in the 3.6/xm channel is about 1.8", the matching of the 
high- 2; galaxy candidates and the IRAC sources was done with a generous, 2" search radius. 
The matched sources were then visually inspected to ensure that the identifications were 
secure. While the rate of the mis-identification depends on a number of factors (such as 
the brightness of both the source and its neighbours, and their proximity), we found that 
roughly one-forth of the identifications in our sample matched unrelated objects. These false 
identifications were thus rejected. Furthermore, we limit our study only to the most reliable 
identifications. In order to avoid ambiguity in interpreting the measured fiuxes, we do not 
include any sources that are blended with foreground objects. Neither do we include any 
sources that are fainter than m^ Q^rn = '26 A mag (the formal S/N = 5 limit for isolated 
point sources in the 3.6/xm channel). In practice, all but two of the matched objects have 
?«3.6/im < 25.3 mag. 

Because of the precise alignment of the IRAC mosaics to the GOODS ACS images, 
the securely matched IRAC sources always have centroids within 0.6" from the centroids as 
measured in the ACS images^. In total, seven z ^ 6 objects (in 3 multiple systems) and 
twenty-two z ^ 5 candidates (3 single objects and 8 multiple systems) were securely matched 



^Scc http: / /data. spitzcr.caltcch.edu/popular/goods/20041027j3nhanccd_vl. 

•^The astrometric accuracy, which is based on brighter point sources, is always good to 0.2 in all channels. 
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with IRAC sources^. Since none of the multiple systems is resolved by the IRAC images, 
from now on we will not distinguish whether an IRAC source is a multiple system or a single 
object, and will simply call it an "object"®. Therefore, our IRAC-detected high-z sample 
consists of three z ^ 6 objects and eleven z ^ 5 objects. Among them, two of the z ^ 6 
galaxies and three of the z ^ 5 galaxies have been spectroscopically confirmed (Dickinson 
et al. 2004; Stanway et al. 2004; Malhotra et al. 2005; Stern et al., in preparation; see 
§3.4). All the objects in our sample are clearly detected in 3.6/im, and some of them are 
detected in 4.5//m as well, where the data are somewhat shallower. The background noise is 
substantially larger and the PSF broader in the 5.8/xm and 8.0//m IRAC channels, and none 
of the galaxies is significant detected at those wavelengths. ACS and IRAC image cutouts 
(and, where available, NICMOS as well) are displayed in Fig. 1 and 2 for the the z 6 and 
z ~ 5 galaxies, respectively. 

Perhaps unsurprisingly, these IRAC-detected objects are also some of the optically 
brightest galaxies in the ACS high-z sample. The faintest Zgso-band magnitude of the three 
z ^ Q sources is 26.88 mag (7^ Tab after adding the fluxes from its two members; detailed in 
§3.2), while the faintest iyrs-band magnitude of the eleven z ^ 5 sources is 26.74 mag (7^48ab 
after combining its two members). Six objects with ^gso < 26.9 in the z ^ 6 sample, and 19 
objects with ^775 < 26.9 in the 2; ~ 5 sample, are not included in the present analysis. Only 
the faintest z !v 6 object and the three faintest z ^ 5 objects are certainly not detected in 
the 3.6//m IRAC data at the 2a level. The others are not included because they are blended 
with other, nearby IRAC sources and were thus rejected during visual inspection regardless 
of whether they were detected. 



3.2. Photometry 

When studying SEDs measured by different instruments at widely separated wavelength 
ranges, one of the largest challenges is to measure a consistent fraction of the light from each 
object in each passband such that the spectral energy distribution is reliable. In our case, 
we utilize matched-aperture photometry as much as possible (using SExtractor of Bertin & 
Arnouts 1996), and choose to bring the measurements as close to the total magnitudes as 
possible. Our current best-effort photometry of these IRAC-detected 2; ~ 6 and 5 galaxy 



''While there seems to be a larger fraction of multiple systems than single objects, we have not yet found 
sufficient evidence that IRAC preferably identified multiple systems. 

^For these multiple systems, we emphasize that they satisfy our high-z color selection criteria either when 
counting their components individually or when combining the components as a whole. 
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candidates is given in Table 1 and 2, respectively. 

For the z ^ 6 candidates, the ACS magnitudes (^775 and 2:850) are taken from YW04, 
which are SExtractor "MAG_AUTO" magnitudes measured through matched apertures us- 
ing the ^850 mosaic as the detection image. The ACS magnitudes {Vqqq, ^775 and Zs5o) of the 
z ^ 5 candidates are also "MAG_AUTO" values, but were obtained using the 2775 mosaic as 
the detection image. For the multiple systems (which are not resolved by either the NICMOS 
images or the IRAC images), their combined magnitudes are listed as well. 

For those sources that are detected in the NICMOS HUDF images, their Jno and 
HiQo magnitudes are based on the NICMOS vs. ACS color indices measured through 
matched apertures using re-binned ACS mosaics as the detection images. The high-resolution 
(0.03"/pixel) ACS mosaics were first block-averaged by 3 x 3 to the resolution of the NIC- 
MOS mosaics (0.09"/pixel), and then registered to the NICMOS mosaics®. SExtractor was 
run using the degraded 2:850 image and ^775 image as the input detection images to extract the 
NICMOS vs. ACS colors for the z 6 and z ^ 5 objects, respectively. These color indices 
were then added to the 2:350 (^775) magnitudes to obtain the final NICMOS magnitudes fisted 
in Table 1 (2). 

Because the angular resolution of the IRAC images is approximately 10-20 x broader 
than that of the HST NICMOS and ACS images, simple rebinning and matched aperture 
photometry is not appficable, and a different method must be used. The GOODS IRAC 
images are very deep, and crowding by neighboring objects can be a problem for photometry, 
even for objects selected to be relatively isolated, like our high- 2; sample. We minimize this by 
using small apertures (3" diameter) and a correction to total flux. The aperture corrections 
were derived from extensive Monte Carlo simulations in which artificial galaxies with varying 
magnitudes and sizes were inserted into the images and then detected and measured using 
the same SExtractor cataloging parameters used for the real data. Because our objects 
are compact, we apply aperture corrections derived for simulated galaxies with half-light 
radii < O'.'S, which average -0.55 mag and -0.60 mag for the 3.6//m and A.B/im channels, 
respectively. The corrected aperture magnitudes for IRAC are given in Tables 1 and 2. 

In addition to the ACS, NICMOS and IRAC data, we also use the deep VLT/ISAAC 
Ks-hand images (Vandame et al., in preparation; Giavalisco et al. 2004a). Only one object, 
#6ab in the z 5 sample, is significantly {S/N > 3) detected. While the inclusion or 
exclusion of this data point is not significant for this object in particular, we use it for the 
sake of completeness. 



^We used the "rotated" (North-up and East-left) NICMOS images included in the version 1.0 public data 
release made by the NICMOS HUDF team on March 9, 2004. 



-8- 



3.3. Assessment of systematic errors in photometry 

The systematic errors in IRAC photometry are largely caused by the source blending. 
To estimate such systematic errors, we used two more methods to estimate the IRAC mag- 
nitudes of these sources. While none of these methods (including the aperture photometry) 
completely eliminates the contamination caused by blending, they are affected differently. 
By comparing these three different kinds of magnitudes, we can estimate the effect of con- 
tamination to these sources. 

One of the additional estimates is the usual "MAG_AUTO" magnitude measured by 
SExtractor. The other one is based on point spread function (PSF) fitting, which is justified 
by the fact that these sources are very nearly point sources at IRAC's resolution. The 
centroids of the sources measured from the ACS images are used as priors to the fitting 
algorithm. The total fitted magnitude of a multiple system is based on the sum of the fitted 
fluxes of the individual components. We use the root-mean-square (rms) dispersion of the 
three types of magnitudes as a measure of the systematic plus random measurement errors 
in these crowded images. The zeropoint calibrations for IRAC are believed to be accurate at 
the 10% level or better, but we allow for additional systematic uncertainty in the absolute 
photometry at this level, and therefore add (by simple summation because it is a systematic 
term) an additional 0.1 mag to obtain the final error estimates for these objects (Tables 1 
and 2). 

For the ACS and the NICMOS photometry, we assess the possible systematic errors by 
comparing our results against those of others. We use the catalogs included in the HUDF 
ACS and NICMOS data release (hereafter "the pubhc catalogs" ) as well as other published 
results in the literature for this comparison. We put the emphasis on the comparison of 
colors because they play a more important role than anything else in the SED analysis. 
Given the different source extraction procedures in the different studies, colors are also the 
most robust metric for such a comparison. 

Since the public ACS catalog is largely an i775-band-based catalog {i.e., the i775-band 

image was used to detect sources and to define source extraction apertures), we only compare 
the photometry of the z ^ 5 sources. Object ^^10 and #14ab are excluded from this 
comparison because they are not in the public ACS catalog. 

For the majority of these sources, both the (Veoe — ^775) and (^775 — ^sso) indices agree 
to within 0.1 mag in these two catalogs. Given that the two catalogs use different fla- 
vors of magnitude (the ACS public catalog uses SExtractor "MAGJSO" while ours uses 
"MAG_AUTO"), and that most of the sources have irregular morphologies, such an agree- 
ment is reasonable, and the size of the discrepancies reflects the impact of choosing different 
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photometric apertures. 

There are three sources whose (Vgoe — ^775) colors differ by > 0.1 mag, namely, #4ab, 
15ab and 48ab. Such discrepancies are mainly caused by the large photometric errors in the 
Veoe-band, and can be explained by the random errors {i.e., those derived based on the S/N 
of the sources) associated with these colors. Only #15ab has (^775 — 2:850) different by > 0.1 
mag. Our measurement gives an index of 0.13 ± 0.04 mag, while the ACS public catalog 
gives 0.28±0.03 mag. This very irregular object has also been studied in detail by Rhoads et 
al. (2005), who used different isophotal apertures on a slightly different version of the HUDF 
images. Those authors give an even more different value of (2775 — Zs5o) = 0.47 ± 0.02 mag 
(after combining their "core" and "plume" components). We believe that the discrepancies 
of different results can also be attributed to the different choices of photometric apertures 
and the contamination from the two physically unrelated close neighbors. 

The means of the color index differences are 0.01 ± 0.12 mag and 0.02 ± 0.08 mag for 
(^606 — '^775) and (^775 — Z85o), respectively. The useful quantities are the dispersions, which 
tell us how these colors can be different because of the systematic errors in their photometry. 
Given the Vgoe-dropout nature of these sources, it is reasonable to assume that the errors 
in both ^775 and 2^850 are comparable, while those in Veoe are larger. Therefore, we assign 
these sources additional errors of 0.056 mag to both 2775 and 2:850, and 0.106 mag to Veoe, 
respectively. These numbers, while added in quadrature, give the corresponding values of 
dispersion in the relevant colors. For the 2; ?a 6 objects, wc assume that their (2775 — 2:850) 
colors suffer a comparable amount of systematic error as the (Vgoe — "^775) colors of the z ^ 5 
sources. Thus we assign them additional errors of 0.106 mag to ^775, and 0.056 mag to 2:850, 
respectively. The errors reported in Table 1 and 2 are the sum of such additional systematic 
error estimates and the random errors based on their S/N. 

We compare our (^775 — 2:359) colors of the three 2: 6 sources to the values (based 
on a fixed-size, 0.5 "-diameter aperture) reported in Bunker et al. (2004), and find that the 
differences can all be explained by the systematic uncertainties that we adopt. 

Since our NICMOS magnitudes are derived based on matched-aperture photometry 
between ACS and NICMOS, we believe those magnitudes are superior to any others for the 
purpose of SED analysis. As these magnitudes are calculated by adding NICMOS vs. ACS 
colors to the ACS magnitudes (i.e., adding Jno — ^775 or Hiqq — ^775 to ^775 for 2; ~ 5 objects, 
and adding Juq — zs5o or Hiqq — zsm to 2:350 for 2; ~ 6 objects), these values should "inherit" 
the systematic errors of the ACS magnitudes. Therefore, we assign an additional error of 
0.056 mag to both bands. Compared to our photometry, the "MAG_AUTO" values in the 
NICMOS public catalog (version 1.0) on average are fainter by 0.12 mag in Jno and 0.24 
mag in i^ieo, respectively. The NICMOS photometry of the three ^ « 6 objects has also been 
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discussed by Stanway et al. (2004), and their values agree to ours to within the quoted errors 
(0.1-0.2 mag). Rhoads et al. (2005) also discussed the NICMOS colors of object #15ab in 
our z ^ 5 sample; while their Jno — Hiqq colors are bluer by 0.19 mag as compared to our 
measurements, such a difference can still be understood if we consider the quoted errors. 

3.4. Galaxies with confirmed redshifts 

Two objects in the z ^ 6 sample have spectroscopic redshift measurements. Object 

#1 was observed with the Keck and VLT observatories and with the ACS grism, and has 
z = 5.829 (Dickinson et al. 2004; Stanway et al. 2004; Vanzella et al. 2005). Based on its 
morphology, the secondary component ^^Ib seems to be a "filament" extending from #la 
(see Fig. 1), and therefore it is taken as at the same redshift as #la. 

The Grism ACS Program for Extragalactic Science (GRAPES; Pirzkal et al. 2004a) 
has confirmed the high- 2; nature of 7^5abc (Malhotra et al. 2005). However, the components 
7^5ab (not resolved by the GRAPES; ID 3450) and ^5c (ID 3503) are identified at rather 
different redshifts of 5.9 and 6.4, respectively. As their contributions to the IRAC fiuxes 
cannot be separated, as a first-order approximation we treat them as at the same, but not- 
yet-determined redshift. We will return to the discussion of the redshift of this system in 
§5.1 again. 

Four objects in the z ^ 5 sample have spectroscopic redshifts. Object #6ab, 7ab and 
15ab have been observed at Keck and are confirmed to be at z = 4.65, 4.78 and 5.49, 
respectively (Stern, in preparation). The GRAPES has also measured z — 4.88 and 5.52 
for #7ab and 15ab, respectively (Xu et al. 2005; see also Rhoads et al. 2005 for an initial 
identification of #15ab, which gave z — 5.42 for this source). Here we adopt the values 
obtained at Keck. The GRAPES has also measured z — 5.05 for object #4ab (Xu et al. 
2005), which we adopt in our discussion. 

4. Overview of Stellar Population Synthesis Models and Procedures 

We use the stellar population synthesis models of Bruzual & Chariot (2003; BC03) to 
analyze the SEDs of these IRAC-detected high-z galaxies. For ease of comparison with other 
studies, we use a Salpeter (1955) initial mass function (IMF). This IMF is incorporated in 
BC03 models with lower and upper mass cut-offs at 0.1 and 100 Mq, respectively. Recent 
work has favored a ffatter IMF below 1 Mq (e.g., Chabricr, 2003). For our purposes, the 
shape of the IMF at low masses only affects the derived stellar mass-to-light ratios (M/ L) 
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in a way that is nearly independent of stellar population age. The Chabrier IMF yields 
masses that are 0.53-0.60 times those derived using a Salpeter IMF; the same will be true 
for galaxies at lower redshifts to which we might compare our results for the high-redshift 
population. Instead, changing the slope or shape of the IMF at high masses can have more 
pronounced effects on galaxy colors and their time evolution, and on the relative mass-to- 
light ratio as a function of stellar population age. E.g., a top-heavy IMF will produce more 
light and have smaller M/L at younger ages, and the luminosity fade more quickly with 
time. As yet, we have little or no observational constraint on the IMF for galaxies at z «i 5 
to 6, and we therefore adopt the Salpeter form to simplify comparison with other published 
work, but systematic uncertainties related to the choice of the IMF should be kept in mind. 
We will return briefly to one IMF-related issue in Sections 5.2 and 6.2. 

The model spectra generated by the BC03 code are shifted across the range of expected 
redshifts with a step size = 0.1, and attenuated by intergalactic HI absorption according 
to the recipe of Madau (1995). The model spectra are then integrated through the bandpass 
response curves and compared to the observations. We consider BC03 models with five 
metallicities Z, from Z0/2OO to solar. We explore a variety of star formation histories (SFHs): 
instantaneous bursts (or simple stellar populations, SSPs); continuous but exponentially 
declining star formation with e-folding timescales from r = 10 Myr to 1 Gyr; and continuous 
star formation at a constant rate. We also consider dust reddening following the empirical 
law of Calzetti et al. (2000) within the range of E{B — V) — 0-2.3 mag. An important 
constraint on the models is that the inferred ages should not exceed the age of the universe 
at the corresponding redshifts. For example, if a set of models give photometric redshift of 
Zp = 6.0 for an object, the legitimate models should only be those that have ages less than 
0.95 Gyr (the age of the universe at 2; = 6 in our adopted cosmology). 

In some cases, the observed photometry cannot be satisfactorily fit by a single-component 
model of the sorts that we have considered. One reason why this may occur is because the 
real SFH is likely more complicated than these basic models. Therefore, we consider the 
superposition of two models when a single model fails to explain the observed SED. This 
"two-component" approach has been used by Yan ct al. (2004) to explain the weak optical 
fiuxes of the IRAC-selected extremely red objects (lEROs). It is similar to the "maximum 
M/L" method of Papovich, Dickinson & Ferguson (2001), where they derived the maximum 
stellar mass allowed by the observed SED of a given galaxy through breaking its light into the 
contributions from a maximally old component (formed at Zf — 00) and from a young com- 
ponent. The differences between the two methods are that introducing one more component 
in our current approach is driven by what is necessary {i.e., when the one-component models 
fail) rather than what is allowed, and that each of the two components can be formed at any 
arbitrary Zf as long as it is larger than the redshift at which a given galaxy is observed. This 
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approach is also similar to that of Berta et al. (2004), but is different in the sense that they 
considered an arbitrary large number of components represent by SSPs, while we limit our 
method to only two components (but not necessarily SSPs) such that the major processes 
can be seen more clearly. 

We carry out the fitting in fiux density units (instead of magnitudes), and minimize 
to find the best-fitting model. When only single-component models are involved, the possible 
free parameters are: redshift (2;), age (T), stellar mass (M), star formation history (r; SSP is 
treated as r = 0), reddening (£'(5 — K)), and metallicity {Z\ When a spectroscopic redshift 
is available, this parameter is fixed. When two-component models have to be evoked, we 
require both components be at the same redshift and have the same metallicity. 

In the practice of minimizing x^, the goodness-of-fit is usually measured by the reduced- 
X^, which is defined as x'^ = x^/ z/, where v is the degrees of freedom of the given problem. In 
our case, v is the number of available bandpasses minus the number of free parameters. The 
objective of our analysis is not only to find the parameters that best describe the SEDs {i.e. 
the ones giving the smallest x'^), but also to find the ranges within which the parameters can 
still give satisfactory fits to the SEDs (ie., the ones giving acceptable x'^). While in principle 
the latter can be achieved by drawing constant x^ boundaries as confidence limits on the 
estimated parameters, in practice it is difficult to take this approach. The fundamental 
difficulty is that the number of available bandpasses is limited, and thus for most of objects 
in our samples the values of v would be zero or negative if we were to vary all the parameters. 

To overcome this difficulty, we take the following approach. For a given object, we 
carry out the fitting procedure in separate runs, and require i/ > 1 in each run. This is 
achieved by "freezing" different parameters in different runs separately. Whenever possible, 
we simultaneously fit for the maximum number of free parameters allowed, i.e., making 
V = \. For example, for an object that has five bands available, we can fit for four free 
parameters at the same time if using single-component models. A possible run then can fit 
for (T, M, z, Z) with (r, E{B - V)) fixed, or (T, M, z, r) with (Z, E{B - V)) fixed, and so 
on. We start from an arbitrary combination of maximum number of free parameters (but 
always including T and M, and also z if it is unknown), and find the best solution in this 
run. Of course, this requires we fix other parameters to some initial values based on our 
best-guess. We deem a fit to be of high quality if it has x'^ <! + (Press et al. 1992). 

When z/ = 1, which is true for most cases, this criterion is x'^ ^ 2.41. If the best-fit in this 
run is a high-quality solution, we then "freeze" the parameters (anything other than T and 
M) to the values of this best-fit, and begin finding the best-fit in a different run involving 
a combination of different parameters. We repeat this process until all free parameters are 
exhausted, and then declare the best solution in the last run as the best-fit model to the 
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object under question. Finally, we alter all the parameters around the best-fit values and see 
to what extent we can still obtain high-quality fits. This approach has a number of caveats, 
one being that such a solution might be a local best-fit but not a global best solution. Here 
we have to neglect such possible caveats. 

5. Stellar populations at z ^ 6 

Here we discuss in detail the fitting results for the three objects in the z ^ 6 sample. 
The passbands that go into the fitting procedure are from ^775 to 3.6/im (see Table 1). The 
best-fit models are shown in Fig. 3, and their parameters are summarized in Tab. 3. 

5.1. Evidence of evolved, massive systems 

We discuss the three z ^ 6 objects in order of easiness in their fitting. 

Object 7^ Tab: This object can be fitted by single-component models. The best-fit 
model (x'^ = 0.66) to this object is a short burst (r = 10 Myr) at Zp = 5.9, with solar 
metaUicity and no reddening. The derived age is 0.1 Gyr and the stellar mass is 4.7 x lO^M©. 
Around this best-fit, we find following constraints on the free parameters: 

(1) r: At < r < 10 Myr (including SSP), we can find high-quality solutions of the same 
quahty (x'^ < 0.7) as the best-fit, and the derived age and mass are also effectively the same. 
The quality of fit is rapidly reduced at r > 10 Myr, although it does not seem to depend on 
T in a linear fashion. We are able to find high-quality solutions sporadically with increasing 
r all the way to r < 70 Myr, and if introducing E[B — V)^ 0.2 mag such solutions can still 
be found until r < 0.1 Gyr. Such solutions with r > 10 Myr, however, seem unstable, as 
changing their ages (T) by only one step in either direction will immediately make x"^ > 5. 

(2) Z: While lowering metaUicities always (although non-linearly) gives at least two- 
times worse x'^) we cannot reject sub-solar metallicties at a significant level. Such models 
generally produce an older age and a larger mass. 

(3) E{B - V): The quahty of fit becomes significantly worse if E{B - V) > 0. At 
E{B — V) — 0.23 mag, any solution has x"^ at least twice as bad as the best fit. If E{B — V) > 
0.45 mag, no high-quality solution can be found. 

(4) T k M: We find high-quality solutions only within 0.05 < T < 0.1 Gyr and 
1.4 X lO'' < M < 1.8 X IO^^Mr. 
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We do not find any acceptable single-component models to fit object #lab and 5abc. 
The best-fit single-component models to #lab and 5abc have x'^ = 10.49 and 2.78, respec- 
tively, and none of them satisfy the high-quality criterion (x'^ < 2.41 in these cases) that we 
adopted. For the sake of completeness, we list the parameters of these fits: (T = 0.9 Gyr, 
r = 0.9 Gyr, M = 1.9 x IO^^Mq, E{B -V) = 0, Z = 0.2Zq) for #lab, and (T = 0.9 Gyr, 
r = 0.5 Gyr, M = 1.6 x IO^^Mq, E{B -V) =0, Z = Q.2Zq) for #5abc. Note that both 
models have a very old age of 0.90 Gyr (the age of the universe at 2; = 5.8 and 5.9 is 0.99 
and 0.97 Gyr, respectively), which means that one has to push the single-component models 
almost to their limits in order to obtain even these mediocre fits. 

Therefore, we have to consider two-component models for these two objects. From 
here on, we use subscripts "E" (evolved) and "Y" (young) to denote the parameters for the 
primary and secondary components, respectively. As we require both components have the 
same redshift and metallicity (see §4), no subscript is put to these two parameters. 

Object #5abc: The best-fit (x'^ = 0.23) to this object has Zp —5.9 and Z — Zq, 
and {Te = 0.9 Gyr, te = 0.2 Gyr, Me = 3.8 x IO^^Mq, E{B - V)e = 0) and {Ty = 0, 
My = 1.0 X 10^Mq,E{B -V)y = 0), respectively. We note that its photometric redshift 
is consistent with the redshift given by the ACS grism (Malhotra et al. 2005; z — 5.9) for 
the component #5ab. Excluding the contribution of #5c (see §3.3) from the SED does not 
change any of our conclusions. 

(I) High-quality solutions can only be found when the young component has Ty < 30 
Myr and E{B — V)y — 0. Its Ty, however, is arbitrary, because any SFH gives effectively 
the same solution at such a young age. For the same reason, this component does not pose 
any constraint on the metallicity. 

(II) If we assume that both components have the same reddening, i.e., E{B — V)e — 0, 
we have following constraints on the primary component: 

(1) T h T: High-quality solutions are available only when te < 0.4 Gyr and Te > 0.4 
Gyr. Thus the age impfies a formation redshift 2;/ > 8.9. The youngest Te (0.4 Gyr) is found 
only when te < 80 Myr. If 0.2 < < 0.4 Gyr, the youngest Te will be 0.7 Gyr, implying 

> 16. 

(2) M: We find high-quality solutions only when 1.3 x 10^° < Me < 6.2 x IO^^Mq. 

(3) Z: At one-step lower metallicity, Z = OAZq, solutions of comparable high quality 
are available. Decreasing metallicity beyond this value steadily decreases the quality of 
fit, and no high-quality solutions can be obtained ii Z — 1/200Zq (the lowest metallicity 
available in BC03). 



- 15 - 



(III) If non-zero E{B — V)e is allowed {i.e., the two components can have different 
reddening), most of the parameters for the primary component are no longer well constrained. 
Even in this case, however, the stellar mass is M > 2.0 x for high-quality solutions. 

We note that none of the solutions with non-zero E{B — V)e has as small x'^ as the best-fit 
model does. If we confine the solutions to solar metallicity, no high-quality solution can be 
found at E{B — V)^; > 1.1 mag. This constraint is largely set by the slope of the SED from 
3.6//m to 4.5/im (larger reddening values would require a much steeper slope). We find that 
most of the high-quality solutions have E{B — V)e'^ 0.45 mag if Z = Zq, which might 
indicate that the true reddening value is indeed insignificant. For low reddening values, the 
high-quality solutions quickly converge to those obtained with E{B — V)e — 0. 

Object #lab (z = 5.83):^ We find no satisfactory one- or two-component model fit 
to this object if we adopt the normal standard of high-quality. As we will discuss in detail 
in the next section, we find that the large x'^ values are mainly caused by its abnormal 
(-2^850 ~ <^iio) color as compared to the BC03 models. This suggests that the models used for 
the secondary (young) component might have some intrinsic shortcomings. As the primary 
(evolved) component is more important to our conclusions, and the models do not have such 
a problem in explaining this component, we still use these models in interpreting this object, 
but loosen our criterion and compare the fits in a relative sense. Therefore, we discuss all 
solutions that have give x"^ < 10- We find a best-fit of x"^ — 5.6, which gives (Tg = 0.5 Gyr, 
te = 0, Me = 3.4 x lO^^M©, E{B - V)e = 0) and (Ty = 0.001 Gyr, My = 2.3 x lO^M©, 
E{B — V)y — 0), respectively. Both components have solar metallicity. 

(I) Solutions of x'^ < 10 can only be found when the young component has Ty < 30 
Myr and E{B — V)y = 0. Similar to the case of #5abc, we have no constraint on the SFH 
of this component. 

(II) If we require the primary component have the same reddening as the secondary one, 
i.e., E{B — V)e = 0, we have following constraints on the primary component: 

(1) T k, T: Solutions of x'^ < 10 can only be found when Te > 0.2 Gyr and te < 0.6 
Gyr. The minimum Te of 0.2 Gyr implies a formation redshift Zf > 7. If 0.5 <te< 0.6 Gyr, 
solutions of comparable quality can only be found at Te > 0.8 Gyr, which means Zf > 20. 

(2) M: Solutions of x'^ < 10 can only be found within 1.4 x 10^° < Me < 5.7 x W^Mq. 

(3) Z: While we see the general trend of worsening fit with decreasing metallicity, 
we are not able to reject models with lower metallicities. However, we point out that if 



^This object is one of the objects discussed by Eylcs ct al. (2005). Our main results broadly agree with 
theirs in the sense that the best-fit stellar mass is > lO^^M© and the age is a few hundred Myrs. 
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Z — 1/200^0, the only x'^ < 10 solutions available are those having very old Tg (> 0.8 
Gyr) and very short te (very close to SSPs). 

(Ill) Allowing E{B — V)e to be different from E{B — V)y makes the models largely 
unconstrained. However, we find that E{B — V)e cannot be larger than 0.4 mag (for the 
same reason as in the case of #5abc) and M cannot be less than 2.0 x 10-'^'^ Mq, or no solution 
of x'^ < 10 can be found. 

The SED analysis for these 2; 6 objects, especially #lab and 5abc, strongly suggests 
that systems as massive as M ~ 10^°Mq already existed when the universe was less than 1 
Gyr old, and they probably started forming their stars several hundred Myr earlier. While 
we do not have a strong constraint on their metallicities, they are consistent with solar 
abundance, implying that they might have already been significantly polluted by metals. 
Dust reddening to the three objects in our sample seems to be moderate or even negligible, 
and the best-fits to their SEDs always have zero reddening. Dust with grayer attenuation 
could, however, cause extinction without reddening. This would also increase the stellar 
masses derived for the galaxies. 



5.2. Very blue rest-frame UV colors and hint of very nicissive steirs 

As mentioned above, ^^lab has an abnormally blue (-2850~<^iio) color as compared to the 
BC03 models^. Such a very blue rest-frame UV color is not unique to this object, however. 
In fact, most of the z f« 6 candidates in YW04 that have NICMOS photometry show such a 
trend in their (^§50 — J\\q) colors, some of which even have (zgso — Jno) ~ —0.5 mag (see Fig. 
1 of YW04; see also Stanway et al. 2004). The matched-aperture photometry described in 
§3.3 further confirms this result, and gives more accurate measurements of the colors. No 
matter what the specific SFH is, a galaxy at this redshift should have a rather fiat SED (in 
terms of ]y vs. A) in the rest-frame UV, i.e., its (-Zgso — ^110) colors should be close to zero. 
The NICMOS Jno bandpass heavily overlaps that of the ACS ^350 bandpass, making it quite 
difficult to produce colors as blue as those observed here. 

Fig. 4 compares the measured (2:850 — J^iio) colors of the z ~ 6 galaxies that have 
NICMOS photometry against the expected (-2850 ~ ^110) colors as a function of redshift as 
derived from a series of solar metallicity models. The models arc drawn from those mentioned 
in §4, and are of four types of representative SFH: SSP, r = 0.1 Gyr, r = 1.0 Gyr, and 



*The discrepancy is not likely caused by a zeropoint error in Jno, as the zeropoints in the NICMOS 
passbands are accurate to 0.01-0.02 mag (R. Thompson, private communication). 
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constant star-formation. The ages of the models are 0.001 Gyr, 0.01 to 0.09 Gyr (0.01 Gyr 

step-size), and 0.1 to 1.1 Gyr (0.1 Gyr step-size). The youngest track is at the bottom, and 
the oldest track is at the top. The three IRAC-detected sources are shown as filled squares: 
#lab is put at z = 5.83 without a redshift error bar, ^^Sabc is put at z = 6.0 with an uneven 
error bar from z = 5.9 to 6.4, and 7^7ab is put at z = 6 with an error bar indicating its 
possible redshift range of 5.5 < 2; < 6.5. The open squares (all put at z = 6 and with error 
bars from z = 5.5 to 6.5) are other z ^ 6 objects that have NICMOS photometry. 

While object ^7ah can be explained reasonably well by a population with age of ~ 
0.1 Gyr, the majority of the sources arc problematic. Regardless of the SFH, the bluest 
color "boundary" is always defined by the zero-age track (almost indistinguishable from the 
O.OOlGyr track in the figure), which cannot extend to (^§50 — ^110) < —0.2 mag at any 
redshift. We emphasize that object #lab, whose redshift is unambiguously known, is not 
consistent with any of these models. Such blue colors cannot be attributed to low metal 
abundance, as the most metal-poor (1/200 of solar, or Z = 0.0001) models of BC03 have 
essentially the same bluest color boundary as shown in Fig. 4. A less severe intervening 
HI absorption than that of Madau (1995) might possibly explain the bluer colors at z > 6 
if (^850 — <^iio) > mag, but cannot explain the {zq^q — Jno) < —0.2 mag colors at z < 6. 
Another possibility would be the presence of a strong Lya emission line in the 2:850 band, 
which could make the (2:850 — ^110) color bluer than the models predict. In order to reproduce 
the observed colors of #lab, however, the Lya line would have to have rest-frame equivalent 
width (EW) > SOOA. This is not the case for #lab; we examined its spectrum obtained at 
the Keck (Dickinson et al. 2004), and found that it only had a moderately strong Lya line 
with rest-frame EW=18A. Neither is it the case for the major component of #5abc, which 
do not show Lya emission line based on the grism spectra of the GRAPES program. 

These very blue rest-frame UV colors might suggest that our model templates do not 
sufficiently populate the entire parameter space. One limitation of the BC03 models is that 
the high-mass end of the adopted Salpeter IMF cuts off at 100 Mq. If more massive stars 
are included, such very blue colors could be explained. For example, a toy galaxy model 
consists of only the hottest stars {T^ff ~ 30, OOOi^T; solar metallicity) from Lejeune et al. 
(1997) can easily explain the observed (2:850 — ^110) colors. The lifetime of such very hot 
(thus very massive) stars, however, is very short (in any case < IMyr), which indicates that 
we might be watching a special episode when these early galaxies were actively forming the 
most massive stars. 

As varying the IMF of the BC03 models is beyond the scope of this paper, here we still 
adopt these models as they are, but caution that the SED fitting based on these models 
will likely have problems in explaining the rest-frame UV part of some of the SEDs. In 
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spite of this drawback, one conclusion seems inescapable from this analysis is that most 
galaxies shown in Fig. 3 (including #lab and #5abc) have young components indicative 
of very recent (< 30 Myr) star formation. Furthermore, the dust reddening in these young 
components is likely minimal, otherwise their intrinsic rest-frame UV colors would be even 
bluer and more difficult to explain. 

6. Stellar populations at z 5 

In this section, wc discuss in detail the fitting results for the z ^ 5 objects. The 
passbands that go into the fitting process are from Vqqq to 4.5/im. The best-fit models are 
shown in Fig. 5 and 7, and their parameters are summarized in Tab. 4. 

6.1. SED fitting results 

We first consider the five objects that have NICMOS photometry. 

Object #6ab {z — 4.65): This object can be well explained by single-component 
models. As its redshift is known, with its 8 bands of photometry the degrees of freedom are 
u — 3. Thus our criterion of high-quahty is x'^ < 1-82. While high-quahty solutions are 
available at all metallicities, the goodness-of-fit improves with smaller Z values. The best-fit 

has x'^ = 0.45, and has (T = 1.3 Gyr, r = 0.4 Gyr, M = 2.8 x IO^^Mq, E{B -V) = 0, Z = 
1/200Zq). Wc find high-quality solutions only within the following ranges: E{B ~V) < 0.23 
mag, T > 0.2 Gyr, T > 0.5 Gyr, and 1.8 x 10^° < M < 5.1 x IO^^Mq. The lower bound of 
its age implies Zf > 7.3. 

Object #15ab (z = 5.49): Our models (either single-component or two-component) 
cannot produce any fit with x"^ < 40 for this object because of its abnormally blue (^759— 2:850) 
color. As discussed in previous section, such an abnormal color of this object in particular 
could be explained by a less severe HI absorption in the IGM along the sight-line. If we ex- 
clude the ^775 band from the fitting process based on the argument that our knowledge about 
the IGM HI absorption is imperfect, we can find high-quality fits using single-component 
models. The best-fit has x'^ = 1-63, and gives (T = 1.0 Gyr, r = 0.6 Gyr, M = 2.2 x 10^°Mq, 
E{B -V) ^0, Z ^ Zq). We find no high-quahty solutions unless E{B -V) ^0. On the 
other hand, we have no obvious constraint on its metallicity. High-quality solutions are 
available if 0.3 < r < 0.7 Gyr, 0.8 < T < 1.0 Gyr, and 1.8 x 10^° < M < 3.6 x IO^^Mq. 
The age range implies that the formation redshift of this object is likely Zf > 15.8. Rhoads 
et al. (2005) considered whether this object could be powered by an active galactic nucleus 
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(AGN) , but did not reach a conclusive answer. We do not find any compelling evidence in 
its SED that this object harbors an AGN. 

Object #48ab: Regardless of the values of other parameters, the photometric redshift 
of this object is well constrained at Zp — 5.2. The best single-component fit to this object 
gives x'^ = 2.48, which does not strictly satisfy our high-quahty criterion (x'^ = 2.41 in this 

case). The parameters of this fit are (T = 0.03 Gyr, r = 0, M = 9.4 x IO^Mq, Z = 1/5^0, 
E{B — V)e = 0). The large x^ of this best-fit comes mostly from the discrepancy in the 
IRAC 3.6/xm band, where the model is lower than the observation by a factor of six. 

Thus we also explore two-component models. We find the following: 

(1) High-quahty fits can be obtained only when E[B — V)y — 0. If we hmit the reddening 
of the primary component also to E{B — V)e — ^, the fit only improves moderately. The 
best-fit (x'^ = 1.82) has a metallicity of Z = 0.4^0, and has {Te = 1.0 Gyr, te = 0, 
Me = 1.2 X IO^^Mq) and (TV = 10 Myr, ry = 0, My = 2.2 x 1O^M0). 

(2) If we allow non-zero E{B — V)e, very high-quality fits (x'^ ^ 0.3) can be obtained. 
However, parameters of these fits span vastly different ranges and thus are not constrained. 
For example, the fitted stellar mass Me can be arbitrarily high with increasingly severe 
reddening. Its age is also arbitrary with different reddening values, spanning from 10 Myr 
to 1.0 Gyr. 

Therefore, we consider the fit to this object uncertain at this point. 

The best-fit models of the above three objects are shown in Fig. 5. The remaining 
two objects that have NICMOS photometry, #4ab [z — 5.05) and Tab {z — 4.78), cannot 
be fitted by any models (single- or two-component) that we considered. The discrepancy 
between the models and the observations is mainly in the rest-frame UV, which is very likely 
a problem intrinsic to the current BC03 models (for example, the cut-off at the high-mass 
end of the adopted IMF) rather than a problem that can be attributed to more complicated 
SFH. 

Objects without NICMOS photometry can all be well-fitted by single-component mod- 
els. As we do not have Jno and i^ieo bands that arc critical in determining the rest-frame 
UV slope and UV-to-optical slope, the parameters (except their photometric redshifts, which 
are largely determined by their Vqqq — ij-j^ colors) are not tightly constrained for most ob- 
jects. The details of the fits to these objects are given below. Object ^^18 is completely 
unconstrained the model fitting. It has comparatively red colors in both Veoe — "^775 and 
^775 — z^bo-i suggesting that it may lie near the upper end of the redshift range for the z ~ 5 
color selection. With a very blue (zgso — m(3.6//m)) = 0.13, it is difficult to further constrain 
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the stellar population parameters of the object in the absence of additional near-infrared 
photometry. While this object is not further discussed here, its SED is shown in Fig. 8 for 
completeness. 

Object #10: The best-fit (x^ = 0.51) to this object has {zp = 5.2, T = 0.4 Gyr, 
r = 0.2 Gyr, M = 9.8 x 10^ Mq, E{B -V) ^ 0, Z ^ Zq). High-quality solutions can be 
found if {E(B -V) < 0.23 mag, T > 20 Myr, 1.1 x 10*^ < M < 4.8 x IO^^Mq). However, we 
do not have any constraint on either Z or r. 

Object #13: This object has a best fit of = 1.0 with parameters of {zp — 4.7, 
T = 1.3 Gyr, r = 0.4 Gyr, M = 1.8 x IO^^Mq, E{B -V)^0, Z^ I/2OOZ0). High-quality 
solutions are available only when {Z < Zq, E[B — V) < 0.23 mag, T > 0.5 Gyr, r > 0.2Gyr, 
0.8 X 10^ < M < 3.8 X 101°Mq). 

Object #14ab: The best-fit to this object has = 0.11, and its parameters are 
{zp = 4.5, T = 0.7 Gyr, r = 0.2 Gyr, M = 2.1 x W^Mq, E{B - V) = 0, Z = Zq). 
High-quality solutions can only be found when {E{B — V) < 0.23 mag, 0.09 < r < 0.8Gyr, 
T > 80 Myr, 1.3 x 10^° < M < 7.8 x 10^° Mq). However, we do not have constraint on its 
metallicity. 

Object #20: This object can also be very well fitted. The best-fit has = 0.14 with 
parameters of (zp = 4.9, T = 0.9 Gyr, r = 0.7 Gyr, M = 3.7 x IO'^Mq, E{B ~ V) = 0, 
Z = I/2OOZ0). However, we do not have any constraint on either T or r, as high-quality 
solutions can be found at any vahics. The metallicity is also unconstrained, but there is 
evidence that the goodness-of-fit improves towards lower Z. The estimate on stellar mass 
spans a wide range, from 5.0 x 10^ to 1.0 x lO^^M©. The only tightly constrained parameter is 
the amount of reddening, for which we find E{B — V) < 0.23 mag in order to find high-quality 
solutions. 

Object #51abcde: The best- fit to this rich group has following parameters: {zp = 4.6, 
T = 1.0 Gyr, r = 0.4 Gyr, M = 1.0 x IO^^Mq, E{B - F) = 0, Z = l/200Zo). We do 
not have constraint on either r or Z, but find that the goodness-of-fit improves towards 
lower Z. We can find high-quahty solutions only if {E{B — V) < 0.23 mag, T > 50 Myr, 
2.0 X 10*^ < M < 2.2 X 1O^°M0). 

The analysis above also shows that galaxies with stellar masses of order M ~ lO^^M© 
existed at 2; ~ 5, and that they likely formed at much earlier epochs. Similar to the 2; ~ 6 
objects, the reddening in these z ^ 5 galaxies is moderate or even negligible. While the 
fits to the SEDs of a number of galaxies show a trend of improved quality towards lower 
metallicities (counter-intuitively as compared to the results in the ^ w 6 sample), they are 
all consistent with solar abundance. 
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6.2. Blue rest-frame UV colors in the z 5 sample 

Similar to the case of the blue rest-frame UV colors in the z ^ 6 sample discussed in 
§5.1, we find that all the z ^ 5 objects have blue (2775 — -2350) colors indicative of recent 
star- formation, and that some of them also have abnormally blue (^775 — -2350) colors that are 
hard to explain even the zero-age models of BC03. This is illustrated in Fig. 6, where the 
observed data points are superposed on the "isochrone" tracks like those shown in Fig. 3. 
The four objects with confirmed redshifts are the squares without horizontal error bars, while 
the remaining objects have redshift error bars indicating a possible range of 4.5 < z < 5.5. 

The objects that have the largest deviations from the models are #4ab and #15ab, 
which are at 2; = 5.01 and 5.49, respectively. The discrepancies, however, are different for 
these two objects, and might represent two different mechanisms. Object #4ab deviates 
from the models similarly to the object ^lab in the z ^ 6 sample, and such a discrepancy 
might be due to the high-mass end cut-off of the IMF adopted by BC03. The blue color of 
object 7^15ab, on the other hand, cannot be explained in this way. Although it might be 
reproduced if the galaxy redshift were significantly smaller, we note (§3.4) that the redshifts 
from both Keck and the GRAPES ACS grism spectrum agree very well {z = 5.40 and 5.52, 
respectively), so we consider this explanation unlikely. Given its redshift, its intrinsic UV 
emission in the iy/s-band, no matter how strong it is, would be easily quenched by the 
LycK absorptions of the IGM HI clouds along the sight-line. The fact that its observed ^775- 
band flux is much stronger than expected suggests that the IGM HI absorption along this 
particular sight-line might not be as strong as what Madau (1995) prescribed. If we keep the 
absorption scheme of Madau (1995), the photometric redshift for this object would be well 
constrained at Zp = 4.60 instead of the spectroscopic redshift z = 5.49. A toy model that 
offers the simplest explanation to its observed color, therefore, is that the distribution of the 
IGM HI clouds still follows the law of Madau (1995) all the way to z !v 4.6, but there are 
very few H I clouds from z ^ 4.6 to 5.49. We emphasize that such a toy model might not be 
physical, and that a plausible model should wait until sufficiently high S/N, high resolution 
spectrum is available for this object. 



7. The number density of massive galeixies at ^ 5 to 6 

The analysis presented here shows that galaxies with stellar masses M > lO^^M© existed 
at 2; f=:i 5 to 6, when the universe was only ~ 1 Gyr old. Moreover, some of those galaxies 
appear to have been forming stars for several hundred Myr, beginning at quite large redshifts, 
z ^ 7 to 20. Eyles et al. (2005) derived similar properties for the two z 6 galaxies that they 
analyzed. The stellar masses of these galaxies are similar to those that have been estimated 
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for typical ("L*") Lyman break galaxies ai z ^ 3 (Papovich et al. 2001), when the universe 
was roughly twice as old as it was at 2; ~ 5-6. As we have noted in §3.1, the high-z galaxies 
that are detected by IRAC are also among the optically brightest star-forming galaxies at 
these redshifts, and thus may be at the upper end of the mass function, although this has 
yet to be properly quantified. We may ask, however, if the existence of objects with these 
masses is consistent with predictions from models of galaxy formation in a A-dominated 
cold dark matter (ACDM) universe. There has been a lively debate in the literature as 
to whether such models can reproduce the observed number density and stellar population 
properties of massive galaxies at lower redshifts {z ^ 1 to 3, e.g., Fontana et al. 2003; Poll 
et al. 2003; Nagamine et al. 2004; Somerville et al. 2004). While our HUDF+IRAC sample 
is not complete in any sense, it is still useful to compare the lower limits on massive galaxy 
space density thus derived against the models in order to extend these tests to earlier cosmic 
times. 

For this purpose, we compare our data to recent hydrodynamic simulations in a ACDM 
universe as described in Nagamine et al. (2004) and Night et al. (2005). Those authors kindly 
provided us with mass functions from two types of simulations, namely, a Smoothed Particle 
Hydrodynamics (SPH) simulation and a Total Variation Diminishing (TVD) simulation. 
The SPH model simulates a co-moving volume nearly 100 times larger than that of the TVD 
model (box sizes of 100h~^ Mpc and 22h~^ Mpc, respectively), and therefore provides better 
statistics for rarer, high mass objects. In the mass range of overlap, however, and before 
resolution effects limit the SPH simulation, the TVD model predicts roughly two to three 
times as many galaxies per unit volume at fixed stellar mass. 

Our galaxy sample is by no means complete due to the fact that we have considered only 
objects that were reasonably isolated in the IRAC images and therefore whose detection and 
photometry was not subject to blending issues. We may nevertheless estimate the stellar 
mass limit at which our sample should be complete, modulo IRAC crowding. We have not 
rigorously quantified the effects of incompleteness. However, in §3.1 we have noted that this 
relatively isolated, IRAC-detected subsample makes up about one third of the total number 
of galaxies in each color-selected sample range down to the ACS magnitude limits of 26.9 
in 2^850 (for z ^ 6) or ^775 (for z 5), and that only a few of the other two-thirds fraction 
not analyzed here are clearly undetected in the IRAC data. We infer, therefore, that the 
incompleteness correction in the IRAC-detected sample is unlikely to be more than a factor 
of ~ 3. 

As described in §3.1, object detection in IRAC was done on the weighted sum of the 
3.6//m and 4.5//m images. The 3.6/xm data are deeper and contribute most of to the signal-to- 
noise ratio of this sum, and therefore we consider that band only for the IRAC detection limit. 
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While the formal S/N — 5 point source limit of the GOODS IRAC data is uissum — 26.4 
mag, we will adopt a brighter S/N — 10 limit ma.e^m = 25.65 mag for the present analysis, 
where detection and photometry are more secure, and where the effects of crowding and 
incompleteness are reduced. All but one of the galaxies in the subsample for which we have 
derived mass estimates (#20 at Zp = 4.9) are brighter than this limit. For the redshift 
range covered by the ACS color selection (roughly 4.5 < z < 6.5), the 3.6/im band measures 
rest-frame light in the 0.48-0. 65//m rest frame wavelength range, i.e., on average around the 
rest-frame y-band. Galaxies in the sample must have sufficient ongoing star formation to 
produce the UV hght which allows them to be detected by ACS. Their ages, star formation 
histories, and dust content result in different stellar mass-to-hght ratios (M/L). In practice, 
for the stellar population models which best fit the objects in our sample, we find 0.17 < 
M^/Ly < 0.66 in solar units. The upper bound is not far from the maximum M / Ly that 
we would expect for an unreddened stellar population at this redshift, given the age of the 
universe (1 Gyr aX z = 5.8, yielding {M/Ly)max = 0.77). Therefore, we adopt the upper 
bound of the model-fit M/L range as appropriate for the sort of UV-selected objects in this 
sample, but it would make little difference to use the maximum allowable M/L. At 2; = 6.5, 
for this upper bound to M/ L, our adopted IRAC magnitude limit corresponds to a stellar 
mass limit of 1.6 x lO^^M©. We therefore expect to be able to detect any galaxy above this 
mass limit in the IRAC data at these redshifts, and we adopt this threshold for comparison 
to the models. 

To calculate the effective volumes over which we identify galaxies in our sample, we 
use the sohd angle of the HUDF, and the redshift selection function for our color criteria 
(§2), which we evaluate from Monte Carlo simulations following the procedures described 
in Giavalisco et al. (2004b). We selected galaxies within the clean area of the ACS HUDF 
image that remains after trimming off the noisy edges, and which subtends 10.34 arcmin^. In 
the Monte Carlo simulations, artificial galaxies with realistic distributions of UV rest-frame 
colors, sizes, and magnitudes, and spanning the range of redshifts of interest, were added to 
the real HUDF ACS images and recovered with SExtractor. The fraction of input objects 
that are recovered and which meet the color selection criteria, as a function of redshift, gives 
the redshift selection function, which we then integrate to get the effective volume for the 
sample. Within the HUDF solid angle, our selection criteria correspond to effective survey 
volumes of 2.5 x 10^ Mpc^ and 2.8 x 10^ Mpc^ for the 2; «i 6 and z 5 samples, respectively. 

We use the best-fit stellar masses for the galaxies we have analyzed, as reported in 
Tables 3 and 4. In the 2; 6 sample, there are two galaxies (out of three total) with 
M > 1.6 x IO^^Mq, corresponding to a space density of 0.8 x 10~^ Mpc~^ for this mass 
range. In the z ^ 5 sample, there are four galaxies (out of eight with derived masses) with 
M > 1.6 X IO^^Mq, or a space density of 1.4 x 10~^ Mpc~^. For comparison, at ^ = 6, 
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the SPH and TVD simulations give cumulative number densities of 2.1 x 10~^ Mpc~^ and 
3.7 X 10~^ Mpc""', respectively, or 2.6 and 4.6 times larger than the observed number. At 
z = 5 the number densities in the models are 9.0 x 10"*^ Mpc~'^ and 18.6 x 10~^ Mpc~^, or 6.4 
and 13.3 times larger than the number of galaxies in our sample above this mass threshold. 

Even allowing for a generous degree of incompleteness in our sample, we conclude that 
ACDM models like those considered here can successfully produce enough galaxies at these 
redshifts with masses like those we have inferred. The basic nature of the result is largely 
insensitive to the choices we have adopted such as the IRAC magnitude limit or the maximum 
allowable M/L threshold, unless there were a substantial number of comparably massive 
galaxies with UV fluxes so faint as to be undetected in the ACS HUDF images, or, at least, 
not selected by our color criteria. The actual ratios of the observed to predicted number 
densities are, however, quite sensitive to these details, due to the very steep mass function 
in the simulations {dN/dM^ oc M", with a ^ —2.5 in this mass range, which we note is 
divergent if extended to arbitrarily low masses, below the resolution of the simulations). We 
also note that these models predict quite steep redshift evolution in the number density of 
objects in this mass range (roughly a factor of 4 to 5 increase from ^ = 6 to ^ = 5). Future 
analysis of the GOODS Spitzer data should allow us to test these predictions by examining 
complete, infrared-selected samples that will more nearly approximate populations limited 
by stellar mass at different redshifts. Two ultradeep IRAC fields have also been observed 
in GOODS-N, which will provide additional dynamic range to help test the steepness of the 
mass function predicted by the simulations. Analysis of other data will also help to control 
the potential field-to-field variations in the observed number density that may result from 
clustering, which we may expect would be strong for galaxies at the upper range of the mass 
function at these redshifts. 



8. Summary 

Using the first epoch of GOODS IRAC observations, we have identified mid-infrared 
counterparts of the z ^ 6 and 5 galaxy candidates selected in the HUDF. Six of these 
galaxies have spectroscopic confirmation, verifying that the IRAC instrument can indeed 
probe galaxies to such high redshifts, and as faint as ~ 26.9 mag in the ACS zg^o or iyso 
bands. In this paper, we study three and eleven objects in the z ^ 6 and 5 samples, 
respectively, all of which are reasonably isolated so that crowding in the IRAC images does 
not compromise their photometry. 

Combining the photometry from ACS, NICMOS (in one occasion also ISAAC at the 
VLT) and IRAC, we analyze the rest- frame UV to optical spectral energy distributions of 
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these galaxies and compare them to stellar population synthesis models. IRAC samples 
redshifted rest-frame optical wavelengths where there is a greater contribution of light from 
intermediate-age and long-lived stars, and where the effects of dust extinction are reduced 
compared to those in the ultraviolet rest frame. These data therefore help us to estimate 
the stellar masses and ages of these galaxies. Deep near-infrared data in two passbands from 
NICMOS are also available for many of the objects. These provide measurements of the UV 
luminosity and spectral slope, which are sensitive to the on-going star formation rate and 
degree of dust reddening, quantities that are insufficiently constrained from the ACS optical 
data alone. 

While the SED analysis reveals diversity in the stellar populations of these galaxies, it 
also demonstrates a number of common properties. The following are the most important 
results from this study. 

(1) Galaxies as massive as ~ lO^^M© already existed when the universe was about a 
biUion years old. These stellar masses are similar to those of typical (L*) Lyman- break 
galaxies at 2; ~ 3, when the Universe was roughly twice as old as it was at 2; ~ 5 to 6. 
Two out of the three objects in the z 6 sample, and at least four out eleven objects from 
the z ^ 5 sample, have best-fit stellar masses M > 1.6 x lO^^M©. While the acceptable 
fits to these objects span a considerable range in the stellar population model parameter 
space, the lower bound to their stellar masses it well constrained to be > IO^^Mq. Although 
the galaxy samples arc incomplete, the resulting lower limits on their space density at these 
stellar masses can be comfortably accommodated by at least one set of recent ACDM models 
for galaxy formation (Nagamine et al. 2004; Night et al. 2005). Those models predict a 
very steep mass function and quite rapid redshift evolution of the number density of galaxies 
at fixed stellar mass, predictions which may be tested by the analysis of larger and deeper 
samples from the GOODS Spitzer data set. 

(2) The photometry for most of the galaxies studied here shows evidence for a pro- 
nounced increase in flux density between the rest-frame UV wavelengths (ACS and NIC- 
MOS) and the optical hght measured by IRAC. This is naturally reproduced by a Balmer 
break that results from the dominant presence of stars with ages of a few hundred Myr. 
All three z ^ 6 objects have best-fit ages ^0.1 Gyr, while at least eight z ~ 5 objects 
have best-fit ages > 0.4 Gyr. In particular, the ^ ~ 6 sample strongly indicates that the 
universe was already forming galaxies as massive as ~ IO^^Mq at z > 7, and possibly even 
at z ~ 20. First-year results from WMAP suggest that the reionization of the universe first 
began at z ^ 15-20 (Spergel et al. 2003). The stellar population properties that we derive 
are in qualitative agreement with this picture, demonstrating that stars, which provide an 
important - if not the only - source of reionizing photons, could indeed be formed at such 
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an early epoch. 

(3) While there is no firm constraint on the metallicity of these galaxies, all the objects 
analyzed here are consistent with solar abundance, which indicates that the galaxies in such 
an early stage of the universe might have already been significantly polluted by metals. In 
fact, the best-fit models to all the three z ^ 6 objects have solar metaUicity, and for one 
object there is strong evidence that its metallicity must be higher than 1/200Zq. If these 
metallicities are representative at z ~ 6, this may cast doubt on the suggestion {e.g., Stiavelli 
et al. 2004) that extremely metal-poor galaxies could be the major sources that finished the 
reionization at z ~ 6. 

(4) In all cases, the best-fit stellar population models have no dust reddening, and 
the allowed models generally have low extinction values. Even so, several galaxies, including 
three with spectroscopic confirmation, show abnormally blue rest-frame UV colors compared 
to these unreddened models. We find that in two cases (at z = 5.83 and z = 5.05), an 
IMF extending to very massive stars (> lOOM©), which are not included in the models 
used here, might account for such blue colors. This would have important implications for 
reionization at 2; ~ 6, offering further evidence to justify a large Lyman continuum photon 
escape fraction that has been adopted in a number of studies {e.g., Yan & Windhorst 2004a; 
Yan & Windhorst 2004b). In another case (at z = 5.49), the blue color might be explained 
by a reduced degree of intergalactic HI absorption along the line of sight compared to the 
mean opacity predicted by Madau (1995). 

The authors thanks Dr. Nagamine for providing his most recent numerical simulations. 
We thank the other members of the GOODS team who have contributed to the success of 
the observations and data analysis. We also thank the referee for very helpful comments. 
Support for this work, part of the Spitzer Space Telescope Legacy Science Program, was pro- 
vided by NASA through Contract Number 1224666 issued by the Jet Propulsion Laboratory, 
California Institute of Technology under NASA contract 1407. 
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Fig. 1. — Image stamps of the IRAC-detected z ^ 6 objects. The ACS ^^sso-band "stamps" 
in the left panel are 1.5" x 1.5" in size, and show the morphological details of these objects. 
The images in the other three panels (12" x 12" in size) shows the matching of these sources 
from ACS ^sso to IRAC 3.6fim. NICMOS Juq stamps are also shown when available. The 
id's of these objects are taken from YW04. Object #lab has a spectroscopically measured 
redshift z = 5.83. 
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Fig. 2. — Image stamps of the IRAC-detected 2; ~ 5 objects. The ID's are taken from Yan et 
al. (in preparation). 
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Fig. 2. — Continued. 
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Fig. 2. — Continued. 
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Fig. 2. — Continued. 
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Fig. 3. — The spectra of the best-fit models to the three IRAC-detected z ^ 6 objects are 
shown in this figure. The parameters of these models are also labeled. The observed SEDs are 
plotted as filled squares with error bars. The model spectra are overplotted, and the synthetic 
model photometry integrated through the bandpasses is shown as open squares which often 
overlap with the filled squares. While #7ab can be explained by single-component models, 
both #lab and 5abc require two-component models. For these two objects, the templates of 
the evolved and the young components are shown in red and blue dashed lines, respectively, 
and their corresponding parameters are also labeled in red and blue, respectively. 
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Fig. 4. — Most of the HUDF z ^ 6 objects that have NICMOS photometry show unusually 
blue (^850 ~ <^iio) colors (YW04). The three that have IRAC identifications are shown 
here as filled squares, while the remaining ones are shown as open squares. The one without 
horizontal error bar is object #lab, which has a known redshift of 5.83. The others (including 
#5abc; see explanation in §3.4) are put at z — 6 with error bars indicating their possible 
redshift range of 5.5 < 2; < 6.5. The BC03 models are used to calculate expected colors 
at different redshifts and ages. Four types SFH are considered and are shown separately as 
the "isochrones" in each of the four panels: an instantaneous burst (SSP), continuous but 
declining star formation with e- folding time scales of r = 0.1 and 1 Gyr, and a constant star 
formation. The ages of the models are 0.001 Gyr, 0.01 to 0.09 Gyr (0.01 Gyr step-size), and 
0.1 to 1.1 Gyr (0.1 Gyr step-size). 



- 36 - 



1 r 

0.1 r 

0.01 r 

1 E- 



^0.1 k 



0.01 r 
1 - 



0.1 E- 



0.01 - 



#6ab, z=4.7 
(^'2=0.45) 





J L 



#15ab, z=5.5 
(X'^=1.63*) 



1 



i 



J I L 



#48ab, Zp~5.2 
(X'^=1.82) 



1 



Z=1/200Zq 
T = 0.4Gyr 
E(B-V) = 
T=1.3Gyr 
M = 2.8xlO,>°M 



J I I L 








1 



T = 0.6Gyr 
E(B-V) = 
T=1.0Gyr 
M = 2.2xlO>oM, 



J I I L 




E(B-V),=0 
l.OGyr 

Mv=2.2xl08M^ Me=1.2x101°M. 

J I I L 



0.5 



1 5 
Wavelength (Mm) 



10 



Fig. 5. — The best-fit models to the three z ^ 5 objects that have NICMOS photometry. 
Symbols are the same as in Fig. 3. Object 7^6ab can be well fitted by single-component 
models. Object #15ab, is problematic because of its abnormally blue (^775 — 2:850) color. If 
we ignore i775-band as justified in §6.2 (see also Fig. 6), reasonably high quality fits can 
be obtained using single-component models. Object #48ab, on the other hand, seems need 
two-component models. However, we consider the fit to this object uncertain (see §6.1 for 
details). 
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Fig. 6. — Similar to the abnormally blue rest-frame UV colors observed in some of the z ^ 6 
candidates, a significant number oi z ^ 5 objects in the HUDF also have blue UV colors that 
are inconsistent with any model. For clarity, this figure shows only the z ~ 5 objects in our 
sample. The "isochrone" tracks, as well as the symbols, are the same as in Fig. 4. The four 
objects without horizontal error bars are those that have been spectroscopically confirmed. 
All other objects are plotted at z = 5 with error bars indicating the possible redshift range 
of 4.5 < z< 5.5. 
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Fig. 7. — The observed SEDs of the 2; ^ 5 objects that have no NICMOS photometry are 
shown in this figure, with their best-fit models overplotted. The symbols are the same as in 
Fig. 5. 
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Fig. 8. — Object 7^4ab and 7ab in the z ^ 5 sample do not have good model fit. For 
completeness, their SEDs are shown here. Also shown is #18, which can be fitted by a large 
range of models and thus is unconstrained. 
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Table 1. Photometric properties of the IRAC-detected ^ ~ 6 objects 



ID 


RA & DEC(J2000) 


*775 


2850 


■hit) 


^fl60 






lab 




26.74±0.14 


25.11±0.08 


25.46±0.09 


25.39±0.10 


24.19±0.19 


24.33±0.17 


la 


3:32:40.01 -27:48:15.01 


26.88±0.03 


25.25±0.01 










lb 


3:32:40.04 -27:48:14.54 


29.03±0.18 


27.41±0.07 










5abc 




27.95±0.23 


26.06±0.09 


26.11±0.13 


26.13±0.15 


24.36±0.42 


24.45±0.50 


5a 


3:32:34.29 -27:47:52.80 


29.07±0.21 


26.97±0.05 










5b 


3:32:34.28 -27:47:52.26 


28.98±0.19 


27.17±0.06 










5c 


3:32:34.31 -27:47:53.56 


29.42±0.22 


27.76±0.08 










Tab 




29.62±0.50 


26.88±0.11 


26.43±0.13 


26.36±0.14 


25.30±0.68 




7a 


3:32:37.46 -27:46:32.81 


> 31.0 


27.50±0.07 










7b 


3:32:37.48 -27:46:32.45 


29.70±0.30 


27.78±0.09 











^ The combined magnitudes and photometric errors are given for these objects. These errors include the systematic 
errors (see §3.3). For the ACS passbands, the magnitudes and errors of individual components are also given for reference. 
These errors of individual components in the ACS bands are based on S/N only. 
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Table 2. Photometric properties of the IRAC-detected z 5 objects^ 



ID 


RA & DEC(J2000) 






-^850 






— ma.etim 




4ab 






28.05±0.25 


25.71±0.08 


26.06±0.10 


26.29±0.22 


26.22±0.16 


25.87±0.71 




4a 


3:32:41.09 


-27:46:42.46 


28.14±0.14 


25.75d=0.02 


26.12d=0.04 










4b 


3:32:41.12 


-27:46:43.32 


30.80±0.39 


29.40±0.12 


29.25±0.18 










6ab2 

6a 


3:32:33.98 


-27:48:02.05 


27.53±0.16 

27.63±0.05 


25.81±0.07 

25.88±0.01 


25.69±0.08 
25.73±0.02 


25.66±0.16 


25.47±0.10 


24.11±0.18 


24.35±0.19 


6b 


3:32:34.04 


-27:48:01.84 


30.14±0.29 


28.87±0.10 


29.38±0.28 










7ab 






26.81±0.14 


25.45±0.07 


25.55±0.08 


25.98±0.17 


25.75±0.11 


24.75±0.23 




7a 


3:32:37.96 


-27:47:11.04 


27.31±0.04 


25.90±0.01 


26.03±0.02 










7b 


3:32:37.94 


-27:47:10.82 


27.88±0.07 


26.63±0.02 


26.66±0.04 










10 


3:32:31.38 


-27:48:13.79 


28.61±0.28 


26.23±0.08 


25.71±0.08 






24.72±0.41 




13 


3:32:41.34 


-27:48:43.09 


27.97±0.16 


26.35±0.07 


26.16±0.08 






24.68±0.46 




14ab 






27.51±0.24 


26.17±0.10 


26.09±0.13 






24.01±0.33 


23.93±0.37 


14a 


3:32:37.55 


-27:45:20.59 


27.83±0.14 


26.46±0.05 


26.38±0.07 










14b 


3:32:37.52 


-27:45:20.12 


29.01±0.26 


27.77±0.09 


27.67±0.14 










15ab 






29.06±0.45 


25.86±0.08 


25.73±0.09 


25.74±0.21 


25.82±0.15 


24.69±0.23 


24.21±0.23 


15a 


3:32:33.27 


-27:47:24.94 


29.45±0.34 


26.50±0.03 


26.26±0.04 










15b 


3:32:33.24 


-27:47:24.50 


30.37±0.84 


26.75±0.03 


26.76±0.06 










18ab 






28.14±0.18 


26.10±0.07 


25.21±0.07 






25.08±0.33 


24.51±0.35 


18a 


3:32:40.86 


-27:45:46.19 


28.76±0.09 


26.67±0.02 


25.83±0.01 










18b 


3:32:40.86 


-27:45:46.40 


29.05±0.10 


27.06±0.02 


26.12±0.01 










20 


3:32:45.25 


-27:48:12.46 


28.56±0.21 


26.77±0.08 


26.58±0.09 






26.11±0.50 




48ab 






28.95±0.26 


26.74±0.08 


26.33±0.09 


26.37±0.20 


26.58±0.14 


24.40±0.62 




48a 


3:32:45.80 


-27:47:25.30 


29.71±0.20 


27.34±0.03 


26.88±0.03 










48b 


3:32:45.78 


-27:47:25.26 


29.69±0.23 


27.67±0.04 


27.34±0.05 










51abcde 






27.63±0.17 


26.21±0.08 


26.11±0.09 






25.21±0.37 




51a 


3:32:43.40 


-27:46:26.87 


28.59±0.09 


27.35±0.04 


27.22±0.06 










51b 


3:32:43.36 


-27:46:27.77 


29.50±0.17 


27.74±0.04 


27.68±0.07 










51c 


3:32:43.38 


-27:46:27.80 


29.11±0.11 


27.75±0.04 


27.61±0.06 










51d 


3:32:43.31 


-27:46:28.52 


29.92±0.18 


28.59±0.06 


28.60±0.11 










51e 


3:32:43.44 


-27:46:28.52 


31.49±0.75 


29.46±0.14 


29.28±0.20 











^ Similar to Table 1, but for the IRAC-identified 2 « 5 objects. The object IDs are taken from Yan et al. (2005, in preparation) and 
axe independent from those of the z « 6 objects reported in Table 1. 

^ This objects is detected in Ks band as well: Ks = 24.62 ± 0.24 mag. 
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Table 3. Best-fit results for the IRAC-detected z ^ 6 objects^ 



ID2 


RedshiftS 


TE (Gyr)4 


Te (Gyr)4 


Me (Mq)4 


TY (Myr) 5 


Ty (Myr)5 


My (Mq)5 


Metallicity (Zq) E{B - V) 


lab 


5.83* 





0.5 


3.4x101° 





1.0 


2.3x108 


1 


5abc 


~ 5.9 


0.2 


0.9 


3.8x1010 








1.0x108 


1 


Tab 


~ 5.9 


0.01 


0.1 


4.7x10^ 








1 



1- Parameters of the best-fit stellar population models shown in Fig. 3. 

2- The IDs of these multiple systems are the combination of the IDs of their individual members. 

^ The asterisk indicates spectroscopic redshift, while the leading ~ sign indicates photometric redshift. 
'' These values are for the major, evolved component. 

^' These values are for the secondary, young component. If these values are not present, it means the object can be explained by a 
single-component whose parameters are given to the left. 



Table 4. Best-fit results for the IRAC-detected z ^ 5 objects^ 



ID 


Redshift 


Te (Gyr) 


Te (Gyr) 


Me (Mq) 


Ty (Myr) 


Ty (Myr) 


My {Mq) 


MetaUicity (Zq) 


E{B - V) 


6ab 


4.65* 


0.4 


1.3 


2.8xl0i° 








1/200 





15al) 


5.49* 


0.6 


1.0 


2.2x10"' 








1 





18a,b 


~ 5.2 


I) 


1.0 


1.2xi0i" 





10 


2.2x10'^ 


0.1 





10 


~ 5.2 


0.2 


0.4 


9.8x10^ 








1 





13 


~ 4.7 


0.4 


1.3 


1.8xlOio 








1/200 





14ab 


~ 4.5 


0.2 


0.7 


2.1xl0i° 








1 





20 


~ 4.9 


0.7 


0.9 


3.7x10^ 








1/200 





51abcde 


~4.6 


0.4 


1.0 


1.0x10'° 








1/200 






' Similar to Table 3, but for the ^ ~ 5 sources shown in Fig. 5 and 7. The top part and the bottom part (separated by the horizontal line) 
are for the sources with and without NICMOS measurements, respectively. 



